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Abstract 

Background: Increases in the spatial extent and density of woody plants relative to herbaceous species have been 
observed across many ecosystems. These changes can have large effects on ecosystem carbon stocks and therefore 
are of interest for regional and national carbon inventories and for potential carbon sequestration or management 
activities. However, it is challenging to estimate the effect of woody plant encroachment on carbon because 
aboveground carbon stocks are very heterogeneous spatially and belowground carbon stocks exhibit complex and 
variable responses to changing plant cover. As a result, estimates of carbon stock changes with woody plant cover 
remain highly uncertain. In this study, we use a combination of plot- and remote sensing-based techniques to 
estimate the carbon impacts of pihon and juniper (PJ) encroachment in SE Utah across a variety of spatial scales 
with a specific focus on the role of spatial heterogeneity in carbon estimates. 

Results: At a plot scale (300 m 2 ) areas pinon juniper (PJ) encroached areas had 0.26 kg C rrf 2 less understory 
vegetation carbon compared to un-encroached sites. This lower amount of carbon was offset by an average of 
1.82 kg C rrf 2 higher carbon in PJ vegetation and 0.50 kg rrf 2 of C in PJ surface-litter carbon. Soil mineral carbon 
stocks were unaffected by woody plant cover and density. Aboveground carbon stocks were highly dependent on 
PJ vegetation density. At a 300 m 2 plot-scale, plots with low and high density of PJ forest had 1.40 kg C rrf 2 and 
3.69 kg rrf 2 more carbon than the un-encroached plot. To examine how these 300 m 2 variations influence 
landscape scale C estimates, historical and contemporary aerial photos were analyzed to develop forest density 
maps in order to estimate above ground PJ associated C stock changes in a 25 ha area. This technique yielded an 
average estimate of 1.43 kg rrf 2 of C accumulation with PJ encroachment. Combining this estimate with analysis of 
tree growth increments from dendrochronologies, we estimate that these PJ stands are accumulating aboveground 
C at an annual rate of 0.02 kg C rrf 2 with no slowing of this rate in healthy PJ. This result is in contrast to what has 
been observed in large areas of drought related PJ mortality, where C accumulation has ceased. 

Conclusions: These results illustrate that the encroachment of PJ forests in SE Utah over the last century has 
resulted in a large (and ongoing) accumulation of carbon in PJ trees and surface litter. However, the magnitude of 
the increase depends to on the density of vegetation across the landscape and the health of forest stands. Both 
management activities that remove forest carbon and forest mortality due to drought or wildfire have the potential 
to quickly reverse the multi-decadal accumulation of carbon in these stands. 
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Introduction 

Over the last century, increases in the spatial extent ('ex- 
pansion') and density ('thickening') of woody plants rela- 
tive to herbaceous species has been observed across a 
range of climatic, geologic, and topographic settings [1-3]. 
These physiognomic changes have impacted the structure 
and function of ecosystems in a number of ways. These 
changes include declines in herbaceous species diversity 
and productivity and alterations in the magnitude and 
spatial distribution of soil nutrients and water [3-11]. 
Changes in woody plant extent and density can have a 
large influence on ecosystem carbon storage. In general, 
woody proliferation leads to the accumulation of woody 
carbon biomass and estimates of the carbon uptake asso- 
ciated with woody encroachment suggest these processes 
could be responsible for up to 25% of the total US C sink 
[12-15]. Despite this potential for carbon storage, the im- 
plication of woody biomass changes on carbon stocks are 
highly uncertain and influenced by variety of conditions 
including the vegetation composition prior to encroach- 
ment, soil characteristics, disturbance regimes, land man- 
agement, and climatic factors [16]. As a result there is a 
clear need to better understand the processes associated 
with carbon stock changes during woody encroachment 
and to improve the measurement of carbon stock changes 
in order to better predict how these regional vegetation 
changes affect national C inventories. 

Estimates of changes in ecosystem C with woody en- 
croachment require an understanding of how both above 
and belowground C responds to landscape changes asso- 
ciated with encroachment. However, environmental con- 
trols on above and belowground C may be highly 
variable and, at present, there are no simple relation- 
ships that can be used to predict how ecosystem C will 
change with encroachment [17-19]. Recent studies sug- 
gest that mean annual precipitation (MAP) across a 
range of growth forms (i.e. trees, arborescents, and 
shrubs) is positively correlated with carbon changes fol- 
lowing encroachment. However, this relationship breaks 
down in arid and semi-arid ecosystems where woody en- 
croachment can lead to both sequestration and a net re- 
lease of C [15,17-19]. 

There are two major factors that complicate prediction 
of carbon stock changes following woody encroachment 
in arid and semi-arid ecosystems; first, soil carbon stock 
responses are highly variable across biomes and second, 
the high degree of ecosystem spatial heterogeneity 
within these ecosystems makes it challenging to accur- 
ately estimate large-scale carbon stock changes. Reviews 
of woody encroachment-associated SOC change indicate 
highly variable soil responses ranging from -40 and 
+ 120 percent across a range of US ecoregions [17,19,20]. 
Some of this variation is likely due to climatic and ed- 
aphic factors but some is also almost certainly associated 



with changing spatial patterns of woody and non-woody 
plants in ecosystems undergoing vegetation change. 

Arid and semi-arid ecosystems tend to have heteroge- 
neous and clumped vegetation distributions and the het- 
erogeneity of above and belowground resources, 
including carbon, tend to increase when shrub and for- 
est cover increases in grassland-dominated ecosystems 
[4-7]. Despite this spatial complexity, most studies of 
woody encroachment infer changes in C stocks from 
point measurements of above and below ground re- 
sources that do not incorporate information on the 
spatial distribution, density, or variation in aboveground 
biomass and how these changes relate to SOC 
[4-7,10,21]. Although these point-based measurements 
offer insight into the biogeochemical processes that in- 
fluence carbon stability, they may not adequately capture 
the landscape-scale variation in woody plant density that 
may ultimately control ecosystem C stocks. For the stud- 
ies that have examined the spatial distribution of woody 
plants, patterns of above-ground carbon appear related 
to the size and spatial distribution of tree and shrub can- 
opies whereas belowground carbon can be influenced by 
the presence of canopy and plant stems [21-24]. This 
variability (and associated complexity in sampling 
design) presents a major challenge to ecosystem-scale 
estimates of carbon stock changes following woody 
encroachment. 

In a number of regions across the western US, there 
have been widespread changes in the cover of pinon 
(Pinus edulis, P. monophylla) and juniper {Juniperus 
osteosperma, J. monosperma) (henceforth P-J). Over 
much of the 20th century, these trees have expanded 
into the semi-arid grasslands and shrublands of the 
western US [25-30]. More recently, a multi-year drought 
has lead to widespread (but patchy) mortality of these 
trees, setting the stage for a potential grassland return 
[31-33]. Here we examine the carbon stock changes as- 
sociated with changes in PJ woodland cover in South- 
eastern Utah. We specifically focus on the question of 
how above and belowground carbon stocks change fol- 
lowing PJ encroachment and subsequent drought- 
induced dieback. We then compare plot-level C stock 
measurements to historical and contemporary aerial 
photographs to examine how vegetation heterogeneity 
(density differences) affects estimates of ecosystem scale 
carbon stocks and determine whether object-based 
image segmentation techniques can be used in combin- 
ation with historical photo archives to document re- 
gional changes in woody plant cover and associated 
aboveground C stocks following woody plant expansion. 
Evaluation of belowground C stock changes are only 
conducted on the 0-10 cm depth; previous work investi- 
gating SOC on the sandstone derived soils of study re- 
gion in this research show SOC to be greatest in the top 
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10 cm of the soil profile and that SOC decreases sub- 
stantially with depth [34,35]. Additionally these soils 
have a surface organic layer of little depth (1-3 cm). 

Methods 

Site description 

The study site was located on the Colorado Plateau, spe- 
cifically in the Hop Creek area within the Manti-La Sal 
National Forest on the southwestern flank of Shay 
Mountain (4196245, 622840; UTM NAD 83 zone 12 N). 
The site is located at 2000 m above sea level and charac- 
terized as an Upland Shallow Loam Black Sagebrush 
ecological site [36]. Precipitation in this study site aver- 
aged 393 mm yr" 1 with a mean annual maximum 
temperature of 15°C and mean annual minimum of 2°C 
(http/ /www. wrcc.dri.edu/summary/ climsmut.html) be- 
tween the years 1971-2000. The site has been historic- 
ally grazed by cattle and is currently grazed at an 
intensity of 145 AUM [37]. All analyses for this study 
were conducted on soils derived from the Kayenta Sand- 
stone Formation [38]. Soils were classified as Barx fine 
sandy loam/fine-loamy, mixed, on mesic, Ustollic 
Haplargid [39] (see Table 1 for soil descriptions). 

The Hop Creek area was chosen because comparison of 
historical aerial photographs from 1937 (Figure 1) show 
significant Pinus edulis Engelm. (pinon) and Junipems 
osteosperma(T brr.) (juniper) encroachment into a black 
sagebrush ecological site. The dominant vegetation in the 
unencroached area was comprised primarily of black 
sagebrush and the grasses Achnatherum hymenoides, 
Bouteloua gracilis, and Sporobolus cryptandrus. All field- 
based measurements described below were carried out in 
September 2006 in plots with similar slope position and 
aspect. 

Plot based C stocks 

To quantify changes in C stocks associated with PJ ex- 
pansion into a black sagebrush ecological site we 
established four 300 m 2 plots (15 x 20 m) that differed in 
the degree of tree density and canopy cover. One plot 
was established in an area free of encroachment 
representing the historic sagebrush condition (HSC) 
prior to encroachment. Another plot was established in 
an area with low tree density (LD) and a third plot was 

Table 1 Study area soil characteristics for each study plot 



established in high tree density (HD). Tree density and 
canopy cover differences between the LD and HD plots 
represent the degree to which PJ recruitment and estab- 
lishment occurred over the past century within the Hop 
Creek area (Figure 1). As result, the names LD and HD 
are based on differences in tree density and canopy 
cover. Because we were also interested in the impact of 
the tree mortality on C storage, a fourth plot was 
established where there was significant pinon dieback 
(PJD). 

Within each plot, we estimated cover of sagebrush, 
grasses, litter-duff, bareground, and tree basal canopy 
area (BCA). We also estimated aboveground (sagebrush, 
grass, and tree) and soil (organic and mineral layer) C 
stocks within each plot. Plant cover and above and be- 
lowground C stock changes as a result of PJ encroach- 
ment and PJ mortality were determined by comparing 
the LD, HD, and PJD values to that of the HSC plot. 

To evaluate whether soil C stock spatial distribution 
changes reflect encroaching tree spatial distribution, we 
geo-referenced all trees and soil samples taken from the 
LD, HD and PJD plots and compared their distributions 
to that of the HSC plot. Additionally, trees from each 
plot were cored and dated in order to determine tree C 
accumulation rate within the study plots. Method spe- 
cifics for each measurement are described below. 

Cover estimates 

Fractional cover of sagebrush, grass, litter-duff, and 
bare-ground was estimated visually in 48, 1 m 2 subplots 
within each 300 m plot. These subplots were evenly 
spaced throughout the 15 x 20 m 300 m 2 plot whereby a 
5 m 2 subplot within the entire plot would have four 
1 m subplots for cover estimation. These cover estimate 
subplots were located in the corners of each 5 m 2 sub- 
plot of the 300 m 2 plot. Tree cover estimates were deter- 
mined for the LD, HD, and PJD plots by measuring the 
BCA of every tree in each plot. For large trees, BCA was 
determined by running a field tape around the tree can- 
opy edge to estimate canopy circumference. Canopy cir- 
cumference was then converted to a BCA. For small 
trees BCA was also determined by converting circumfer- 
ence to area. However, small tree circumference was de- 
termined by measuring canopy diameter along two 



Plot 


Sample 
number 


Mineral bulk 
density (g cm" 3 ) 


Duff bulk 
density (g cm" 3 ) 


Duff depth 

(cm) 


Soil C0 3 

(%) 


Clay (%) 


Silt (%) 


Sand (%) 


HSC 


32 


1 .55 (0.02) 


0.1 3 (0.004) 


0.38 (.03) 


1 .97 (0.69) 


7.47 (0.01) 


43.51 (O.09) 


49.01 (0.11) 


LD 


34 


1 .38 (0.02) 


0.25 (0.03) 


1 .53 (0.45) 


1.17 (0.04) 


8.23 (0.08) 


37.35 (0.10) 


54.41 (0.09) 


HD 


41 


1 .29 (0.02) 


0.27 (0.04) 


2.13 (0.35) 


2.99 (0.84) 


9.71 (0.09) 


36.78 (0.11) 


53.01 (0.14) 


PJD 


34 


1 .38 (0.02) 


0.22 (0.02) 


2.09 (046) 


0.96 (0.07) 


8.24 (0.07) 


37.35 (0.11) 


54.40 (0.17) 



All data are calculated means with ± standard errors in brackets. Means for soil texture are based on an N of 6 per plot. Bulk density and soil carbonate values are 
based on the 0-10 cm soil depth. 
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1937 



2006 





Figure 1 Hop Creek study area aerial photos from 1937 and 2006. The 1937 photo shows our study plots. The boxed areas in each photo 
show a 25 ha area where our study plots were located and was the area used for the image segmentation ortho-photo analysis. 



perpendicular lines. All cover estimates are expressed on 
aim 2 basis. 

Carbon measurements 

Carbon stocks of the understory grass and shrub func- 
tional groups was determined by converting cover values 
of each functional group to a standing stock of C. This 
was accomplished in two steps. First, functional group 
specific linear regression equations (Table 2) were devel- 
oped to convert the fractional cover of grass and sage- 
brush into biomass. Linear regression equations used to 
calculate grass (r 2 0.96 p < 0.01) and sagebrush (r 2 0.90 
p < 0.01) biomass were developed by harvesting above- 
ground plant material from eight 1 m 2 plots with a 5 to 
85 percent sagebrush or grass cover range. Harvested 
sagebrush wet weight was recorded in the field for each 
plot. Sub-samples of field weight sagebrush and all grass 
aboveground biomass were transported back to the lab 
and oven dried at 60°C for 48 hours in order to deter- 
mine a dry weight. Dry weight was used to develop the 
linear regression equations for estimating sagebrush and 
grasses biomass as a function of percent cover. One 
square meter plots within 5 m 2 subplots were averaged 
and scaled to represent the cover within the 5 m plot 
that they were located in. These 5 m 2 averages were used 
to calculate sagebrush and grass biomass within each 
300 m plot. Sagebrush and grass biomass was then 
converted into kg C by multiplying biomass by the mean 



percent C content of each plant (Table 2). Average sage- 
brush and grass percent C was determined by homogen- 
izing and grinding five samples of each plant type with a 
shatterbox and analyzing it with an EA 1110 CNS 
combustion analyzer (Thermo Electron Corporation, 
Waltham, MA). Sagebrush and grass kg C estimates for 
each plot were presented on a 1 m 2 basis. 

The methods for determining tree C were based on 
previously published allometric equations [40,41]. To 
calculate tree C, we measured root collar diameter 
(RCD) on all trees within the LD, HD, and PJD plots. 
Root collar diameter is the stem diameter at ground 
level and is used with allometric equations to calculate 
tree biomass [40,41]. For all woodland plots, RCD calcu- 
lated tree biomass was converted to kg C by multiplying 
PJ biomass with plant-specific percent-C content values 
(Table 2). Pinon and juniper C content was determined 
from harvested, oven dried and ground PJ saplings as 
done for black sagebrush and grass. P. edulis and /. 
osteosperma specific allometric models developed by 
Darling in 1967 [40] were used to calculate pinon and 
juniper biomass. For dead trees devoid of needles, the 
green biomass portion of the allometric equation was 
subtracted from dead trees. Green biomass for dead 
trees was estimated using pinon- and juniper-specific 
green biomass allometrics developed by Darling in 1967 
[40]. The Darling 1967 models [40] were used because 
we did not develop RCD tree biomass allometry 



Table 2 Regression equations and r 2 values used to calculate grass, black sagebrush, live PJ and dead PJ biomass and 
the percent carbon (C) within each plant type for determining kilograms of C 



Plant type 




Model 


r 2 


%C 


Grass 


Biomass = 


0.566 * % cover nT 2 


0.96 


40.84 (0.41) 


Black sagebrush 


Biomass = 


12.132 * % cover m~ 2 


0.90 


45.55 (0.51) 


Live Pinon 


Biomass = 


11.41 * (RCD * 0.1 )A 2 - 6664 ) 


NA 


46.14 (0.43) 


Juniper 


Biomass = 


8.256 * (RCD * 0.1 )A 28058 ) 


NA 


45.64 (0.54) 


Dead Pinon 


Biomass = 


(Live Biomass - (1.853 *(RCD * 0.1 )a 20268 )) 


NA 


46.14 (0.43) 



Fernandez et al. Carbon Balance and Management 2013, 8:8 
http://www.cbmjournal.eom/content/8/1/8 



Page 5 of 13 



equations specific to our study area. The Darling allom- 
etry equations were also used because they are for simi- 
lar tree species and were developed in an area similar to 
Hop Creek. We also felt that the Darling equations used 
were more accurate estimators of biomass than those 
developed for estimating biomass at the national scale 
i.e. [42]. 

For soil C estimates, soil samples were collected every 
2.5 m and geo-referenced. However, in plots with PJ en- 
croachment, sampling points were collected 2-2.5 m in 
order to get a better spatial representation of overstory 
C inputs and better assess soil C stocks. The differences 
in sample spacing are addressed in the geostatistical ana- 
lyses described later. At each sampling point, the duff/ 
litter layer was collected when present with a large flat 
spatula and a circular piece of polyvinyl chloride (PVC) 
with a diameter of 10 cm. The depth of each duff/litter 
layer sample was recorded in order to calculate the bulk 
density of samples collected from this layer. Mineral soil 
samples were taken with a volumetric soil sampler (Soil 
Moisture Equipment Soil Core Sampler Model 0200). 
Ten cm long mineral soil samples were collected from 
the point where the duff/litter layer was collected. All 
soil samples were transported back to the lab and oven 
dried at 60°C for 48 hours. After drying, mineral soil 
samples were passed through a 2 mm sieve. This mater- 
ial was then passed through a soil splitter to obtain 2, 
2 g sub-samples, which were then ground to a fine pow- 
der using a mortar and pestle. Sub-samples of soil were 
then acidified with 15% HCL in order to remove carbon- 
ates. Duff samples also had carbonates removed by this 
acidification process. Percent C content for mineral soil 
and duff samples were analyzed with an EA 1110 CNS 
combustion analyzer. All percent C values were 
converted to kg C m' 2 using bulk density values for stat- 
istical and computational analyses. However prior to 
converting to kg C m' all percent C values were carbon- 
ate corrected with soil sample specific carbonate values. 
Soil specific carbonate values were determined by using 
a modified pressure calcimeter method [43]. 

Mineral and duff/litter SOC were analyzed separately 
in order to better understand soil C differences within 
different C pools. Because we wanted to spatially repre- 
sent soil variability in our study, we collected a different 
number of soil samples within each plot and capture 
vegetation heterogeneity, the sample size varied some- 
what across the plots (Table 1). Due to the differences in 
sample size, we used the non-parametric Kruskal-Wallis 
one-way analysis of variance (ANOVA) by ranks to de- 
termine whether SOC differed between plots. The geo- 
metric mean from each Kruskal Wallis test was 
determined to be statistically different by applying the 
multiple mean rank comparison. This analysis should be 
used primarily as a statistical reference for comparison 



of within plot variance vs. the site means. Our design 
does not allow for a true ANOVA analysis (e.g. there is 
only one true plot) and this analysis is used here as if 
our subplot measurements were independent samples 
(noting that many experimental designs use this ap- 
proach without the caveats) and should be interpreted as 
such. All means from the Krustal Walis ANOVAs are 
followed by a bracketed standard error value. 

Spatial distribution of soil C and trees 

The spatial distribution of soil C was evaluated with 
Moran's Index for spatial autocorrelation (MI). Moran's / 
measures the spatial autocorrelation of geo-referenced 
values and evaluates whether their distribution is clus- 
tered, dispersed, or random. The spatial distribution of 
trees was evaluated with the average nearest neighbor 
(ANN) measure. Unlike MI, which measures the spatial 
distribution of values, ANN measures the spatial distri- 
bution of objects. For MI and ANN a value near +1.0 in- 
dicates clustering, a value near -1.0 indicates dispersion, 
and a value near zero indicates a random distribution. 
Each calculated MI and ANN output has a calculated z- 
score. The statistical significance of a clustered or dis- 
persed distribution is determined by comparing z-scores 
to a range of confidence intervals. For example, in order 
for a clustered or dispersed distribution to be statistically 
significant at a p < 0.05 or 0.01, the z-score has to be re- 
spectively between -1.96 and 1.96 or -2.58 and 2.58. 

Tree age 

We evaluated tree age structure of pinons within each of 
the plots to examine tree recruitment rates over the last 
century. Tree recruitment data was used to examine an- 
nual aboveground C accumulation with the occurrence 
of tree encroachment. An increment core was collected 
from all pinons in each plot to determine tree age. We 
did not collect increment cores from /. osteosperma due 
to known difficulties in dating this species. Total tree C 
values were then divided by tree age in order to calculate 
a kg C yr' 1 accumulation rate for each tree. In this ana- 
lysis we use tree data from LD, HD, and PJD plots; how- 
ever, PJD trees were treated as living trees with green 
biomass C at the time of death. Tree cores were col- 
lected from all pinons > 5 cm in diameter. Tree cores 
were collected approximately 10-15 cm above the 
ground. Tree cores were analyzed at the Institute for 
Arctic and Alpine Research's Dendrochronology Lab, 
University of Colorado, Boulder. Increment cores were 
mounted [44] and progressively sanded with FEPA (Fed- 
eration of European Producers of Abrasives) 120-, 220-, 
320- and 400-grit (162, 68, 44.7-47.7 and 33.5-36.5 urn 
respectively [45], sandpaper using a Bosch belt sander 
(Robert Bosch Corp., Farmington Hills, MI.) Tree cores 
were also hand-surfaced using 400 and 1200-grit 
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sandpaper (33.5-36.5 um). Tree cores were then visually 
cross-dated using a previously developed pinon chron- 
ology for NMM (C. Woodhouse, unpublished data). For 
samples difficult to visually crossdate, undated inner sec- 
tions were run against the previously developed chron- 
ology using COFECHA [46]. COFECHA is a software 
program used to check the quality of the cross-dating. 
All matching dates found using the COFECHA program 
were then visually verified. We recorded the first year of 
secondary growth to calculate tree age. Although we 
made 3 to 5 attempts in the field to obtain cores that in- 
cluded pith, we were not always successful. We esti- 
mated dates for cores without pith by overlaying sets of 
concentric circles on the inner rings of the core [47]. To 
correct for the time to coring height, we collected pinon 
seedlings from a range of sites and soil types across the 
region. Average time to a coring height of 15 cm is ap- 
proximately 11 years. Once an inner ring date was 
obtained from the cores, 11 years was subtracted from 
this date to correct for the time to coring height. 

Coupling tree age with aboveground C data, we mod- 
eled annual changes in aboveground C in the LD, HD, 
and PJD plots from the late 1800s through 2005. To esti- 
mate an annual increase in tree C, total tree C was di- 
vided by tree age in order to estimate an average annual 
C increase for each year of growth. Although factors 
such as climate and tree age clearly influence annual 
growth, based on observation of annual growth rings 
there was no simple model (eg. negative exponential) 
that could be applied to annual tree growth at our sites. 
Thus we chose to increase individual tree C by an aver- 
age amount for each year. Since junipers were not dated, 
we used a linear regression equation of tree age as a 
function of root collar diameter from a previous study in 
the region (Barger, unpublished results) to estimate juni- 
per age within each plot (Juniper Age = 2.82(RCD) + 
19.8, r 2 = 0.62). Initial conditions for pre-encroachment 
aboveground C was set at the level of aboveground C es- 
timated from the HSC plots. The starting year for above- 
ground tree C accumulation for each of the plots was 
determined by the oldest tree in each plot (1891, 1906, 
and 1901 for LD, HD, and PJD respectively). In addition 
to accumulation of aboveground C in the tree layer, we 
estimated declines in belowground C in the understory 
layer over time. Tree density in the four plots was a 
strong predictor of understory biomass (Understory 
C = -0.516 * Tree Density + 137.9, r 2 = 0.84). Following 
this, we estimated changes in understory C over time as 
a function of annual historical tree density based on tree 
recruitment into the plots. Historical tree density was 
used to estimate changes in aboveground understory C. 

In this analysis we use tree data from LD, HD, and 
PJD plots. PJD trees were treated as living trees in the 
aboveground biomass accumulation until the year of 



death. After the year of death we assumed that all living 
aboveground biomass in the PJD plots was converted to 
dead biomass. 

Photo analysis of tree C estimate 

The purpose of the aerial photo-analysis was two fold; 1) 
assess the feasibility of utilizing aerial photography for 
the estimating aboveground C over large areas by com- 
paring plot and aerial based estimates and 2) to deter- 
mine tree C accumulation as a result of encroachment 
between 1936 and 2006. To accomplish this task we 
used an object-based image segmentation technique cap- 
able of detecting individual (or clumped) tree canopies 
(henceforth image segmentation) (ENVI Feature Extrac- 
tion module, ITT Visual Information Solutions, Boulder, 
CO) to delineate woody canopies. These canopy areas 
were then imported into a GIS system (ArcGIS v. 9.3, 
ESRI, Redlands, CA) to calculate canopy area (CA) [48]. 
We applied previously published remotely sensed above- 
ground C calculations to calculate tree C from canopy 
area [49]. The spatial extent of this analysis was a 25 ha 
area encompassing our study plots. A single allometric 
tree biomass equation that uses CA and developed by 
Huang et al. [49] was used to calculate tree biomass" 

log 10 biomass(kg) = -4.66 + 1.321og 10 CA(cm 2 ) 

The CA based allometric equation by Huang et al. [49] 
grouped the RCD models developed by Darling in 1964 
[40] for P. edulis and /. osteosperma into one because 
these species have no structural difference when CA is 
used to estimate aboveground biomass in this study re- 
gion. Tree C was calculated by multiplying biomass by 
0.46; the average percent C of pinon and juniper plant 
material. 

To estimate C accumulation in the delineated 25 ha 
area, new canopies were determined by overlaying the 
1937 and 2006 DOQQ aerial photographs. For this ana- 
lysis the quality of 1937 aerial photograph (Figure 1) was 
suitable for visual inspection but insufficient for quanti- 
tative image segmentation analysis. As a result, visible, 
open, woody plant free areas were delineated by on- 
screen digitizing in ArcGIS. The woody plants (and 
aboveground C) derived from 2006 DOQQ located in 
these 1937 woody plant free areas were summarized to 
represent the amount and rate of C accumulated due to 
PJ encroachment within the 25 ha study area. 

The plot-based estimates of aboveground C stocks 
were then compared to the image segmentation esti- 
mates of aboveground C stocks. Comparisons were first 
conducted by extracting image-segmented canopies 
within the LD and HD plots in order to detect tree C 
differences between our plot-based measurements and 
image segmentation derived values. We also scaled plot- 
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based carbon stocks to the 25 ha study area and then 
compared these estimates to the carbon stocks deter- 
mined through image segmentation. 

Results 

Pinon pine was the dominant tree species in the LD and 
HD plots, comprising 88-96% of the tree population 
(Table 3). Pinon age across the study plots supports ob- 
servations from aerial photos that recruitment into the 
open black sagebrush sites occurred in the early 1900s. 
Pinon age data within all plots indicates that trees were 
present as early as 1894 and that a pulse in recruitment 
occurred in the 1940s and 50s with 43-66% of the trees 
across all plots dating to these decades (data not shown). 
Average pinon age (excluding saplings) was similar 
across the PJ encroached plots and ranged from 47 - 
53 years. Tree density in the HD plot was 47 percent 
higher relative to the LD plot (trees ha' 1 , HD = 1567; LD 
= 1067). Tree density in the PJD plot was similar to the 
LD plot with 1000 trees ha' 1 (live + dead) but declined 
to 100 live trees ha' 1 after the regional mortality event in 
the early 2000s (Table 3). Percent tree mortality within 
PJD was 90% and based on outer ring dates all tree mor- 
tality occurred in 2003, with the exception of one tree 
that died in 2004. The HSC plot contained four juniper 
saplings and no pinon trees. 

Evaluation of a 25 ha area in the 1937 orthophoto and 
2006 DOQQ indicates that ~ 32% of the landscape was 
treeless (similar to the HSC). Additionally evaluation of 
the 2006 DOQQ shows that -22% of the encroached 
area had less canopy cover than our HD plot when the 
2006 DOQQ is segregated into 300 m 2 sections 
(Figure 2). We did not attempt to map areas of dead PJ 
woodland in this analysis. 

Vegetation cover and C changes 

Broad changes in understory cover and composition 
were evident with woody encroachment in these plots 
(Table 3). Live tree density was inversely related to 
understory (sagebrush and grass) cover. Sagebrush C 
stocks in the encroached sites were lower and 
corresponded to a decline of 0.1 kg C m' 2 in the LD and 
PJD subplots and a decline of 0.24 kg C m' 2 in the HD 
subplots with the HSC site as reference point (Table 4). 
Grass cover was 11.8 ±0.4% within the HSC sub- 
plots (Table 3) accounting for a C stock of 0.006 ± 
0.002 kg C m' 2 . Grass cover in the LD was similar to 
that of HSC, lower in the HD, and greater the PJD sub- 
plots (Table 4). 

Total aboveground carbon was 140-712% higher n PJ 
encroached subplots relative to HSC subplots. Com- 
pared to the 0.42 kg m' 2 of C in the HSC plot, the LD, 
HD, and PJD plots contained approximately 1.02, 2.99, 
and 0.99 kg more aboveground C per m 2 respectively. 



These changes occurred despite declines in grass and 
sagebrush C. 

Litter and soil C stocks 

Total soil C stocks were 32-34% higher in encroached 
subplots, which was driven by differences in the amount 
of C in the duff/litter layer (Table 4). The difference in 
duff/litter C is attributed to PJ encroachment because 
duff/litter in encroached subplots is comprised primarily 
of PJ needles and small PJ twigs. There were significant 
differences in duff/litter C between subplots (Kruskal 
Wallis: H = 15.56, p = 0.001). A multiple comparison of 
mean ranks for all groups test showed duff/litter C on a 
kg m' 2 to be lower in the HSC subplots than the LD 
subplots and that duff/litter C is highest in the HD and 
PJD subplots. 

Carbon stocks in the surface 0-10 cm of the mineral 
soil layer followed a different pattern from the duff/litter 
soil layer. Differences in mineral SOC did exist between 
the subplots but not in a consistent pattern and the 
magnitude of the differences were relatively small 
(Table 4; Kruskal Wallis: H = 12.68, df=3, n=141, p = 
0.005). Mineral SOC respectively averaged 1.52 ± 0.06, 
1.68 ±0.07, 1.37 ±0.05, and 1.44 ±0.06 kg m' 2 in the 
HSC, LD, HD, and PJD subplots. 

Spatial C distribution patterns 

Pinon and juniper encroachment influenced the spatial 
distribution of the duff/litter but not the distribution of 
mineral layer SOC. Mineral layer SOC was randomly 
distributed in all plots (Table 5) regardless of above- 
ground vegetation cover or condition. Duff/litter C dis- 
tribution was determined to be random in the HSC plot 
but clustered in the LD, HD and PJD plots. The clus- 
tered duff/litter distribution occurred in association to 
the spatial distribution of trees. 

Aboveground C changes over time 

Estimates of aboveground C stock changes in the 
encroached plots over the last century suggests that 
these upland PJ sites have a high capacity to rapidly ac- 
cumulate aboveground C. Modeled aboveground C ac- 
cumulation based on tree recruitment data suggests that 
aboveground C has exponentially increased over the past 
century (Figure 3). These exponential increases in above- 
ground C, however, did vary between the LD and PJD 
subplots relative to the HD subplots; a pattern likely 
driven by differences in tree density. 

Image segmentation photo-analysis tree carbon 

Comparison of the 1937 photo and 2006 DOQQ shows 
that in 2006 PJ occupied ~68% of area that was open in 
1937 (Figure 1). Image segmentation analysis delineated 
nearly the same number of canopies and canopy area as 



Table 3 Total number of trees per plot broken down by species, mean tree age, mean basal canopy area (BCA), tree trunk basal diameter (TCD), and percent 
grass, sage, duff/litter, and bareground cover per m 2 



Plot 


Live tree density 
(Stems ha" 1 ) above 
Number of trees below 


Mean tree 
Age (yrs) 


Mean tree BCA 

(m 2 ) 


Mean TCD 
(cm) 


Mean grass 
cover (%) 


Mean sagebrush 
cover (%) 


Mean duff/litter 
cover (%) 


Mean bareground 
cover (%) 


HSC 


133 4 (J saplings) 


NA 


NA 


NA 


1 1 .85 (043) 


37.03 (0.95) 


NA 


41.38 


LD 


1033 31 (30 P, 1 J) 


51 (3) 


4.2 (0.7) 


10.64 (1.04) 


1 1 .97 (243) 


25.62 (2.57) 


27.60 (4.06) 


31.98 (4.04) 


HD 


1400 42 (37 P, 5 J) 


53 (2) 


5.1(0.8) 


13.70 (1.32) 


2.71 (1.06) 


14.16 (1.88) 


53.84 (5.47) 


27.08 (4.42) 


PJD 


100 30 (15 P, 3 J) 


47 (1) 


4.9 (2.5) 


8.54 (1.55) 


20.41 (2.32) 


25.83 (2.25) 


25.21(3.36) 


24.58 (3.22) 



For the HSC plot duff/litter and bare ground estimates are combined and presented in the bareground column. Percent cover estimates do not add up to 100 because forbs and shrubs other than black sagebrush are 
not included. Values in brackets are ± 1 standard error of the mean. P = Pinus edulis, J = Juniperus osteosperma. 
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Figure 2 2006 1-m resolution DOQQ of the of the Hop Creek area segmented into 300 m cells showing canopy area variance. The 

color ramp represents m 2 canopy cover for grid cells lying within the 25 ha image segmented area. 



measured on the ground in the LD and HD plots 
(Table 6). This image segmented canopy area when 
converted to a C stock with canopy based allometry re- 
sults in a stock estimate that is approximately 30% 
greater than what was directly measured in the LD plot 
and a 10% lower in the HD plot. 

Comparison of ground and image-segmentation-based 
estimates of C/tree yielded very similar estimates. Tree 
C within the LD and HD plots measured 17.2 kg with 
ground based measurements and 17.7 kg with image 
segmentation. The differences between the methods are 
evident when scaling the directly measured values from 
the 300 m field plots to the larger 25 ha study area. A 
direct scaling of the LD and HD plots to 25 ha yielded 
estimates of 267 and 781 Mg C respectively. In contrast 
image segmentation analysis, which measures the entire 
25 ha area yielded an estimate of 358.93 Mg C. 



Differences in the plot-level and remote sensing mea- 
surements occur because plot measures cannot account 
for the tree density variation in the larger woodland area; 
whereas the image segmented C estimate does incorpor- 
ate canopy density woodland variation (Figure 2). 

Combining the remote sensing estimate C stocks with 
the tree ring analysis, we can estimate the rate of C in- 
crease in this area. We estimate that 869 trees ha" 1 ex- 
panded into treeless areas between 1937 and 2006 
(Figure 4) resulting in a woody C increase of 15.36 
Mg C ha' 1 ; and annual C accumulation increase of 
0.22 Mg ha' 1 . 

Discussion 

The expansion of pinon-juniper populations into ecosys- 
tems in SE Utah, USA, causes a predictable decline in 
grass and shrub cover in proportion to tree canopy 



Table 4 Total and individual kg m 2 C stocks for each plot 


Carbon Stock 


HSC (kg m" 2 ) 


LD (kg rrf 2 ) 


HD (kg m" 2 ) 


PJD (kg m 2 ) 


Pihon and Juniper 


0.0003 (NA) 


1.13 (NA) 


3.24 (NA) 


1.10 (NA) 


Black sagebrush 


0.41 (.01) 


0.31 (0.03) 


0.17 (0.02) 


0.31 (0.02) 


Grass 


0.006 (.002) 


0.006 (0.001) 


0.001 (0.0006) 


0.011 (0.001) 


Duff/Litter 


0.05 (0.01) c 


0.42 (0.11) b 


0.71 (0.13) a 


0.53 (0.14) a 


Mineral Soil 


1 .52 (0.06) ab 


1 .68 (0.07) a 


1.37 (0.05) b 


1 .44 (0.06) b 


Total C 


1.99 


3.55 


5.49 


3.39 



Numbers in brackets are ± standard error values. For duff/litter and mineral soil values with the same letter are not significantly different. Pihon and Juniper C 
stocks are calculated with the allometric equations developed by Darling [40]. 



Fernandez ef al. Carbon Balance and Management 2013, 8:8 
http://www.cbmjournal.eom/content/8/1/8 



Page 1 0 of 1 3 



Table 5 The spatial distribution of mineral and the duff/ 
litter soil carbon (C) and encroached plot trees 



Plot 


HSC 


LD 


HD 


PJD 


Ml mineral soil C 


-0.15 


-0.01 


-0.13 


-0.01 


Z-score Ml mineral soil C 


-0.88 


0.09 


-0.86 


0.38 


P-value Ml mineral soil C 


NA 


NA 


NA 


NA 


Distribution mineral soil C 


random 


random 


random 


random 


Ml duff/litter C 


-0.04 


0.32 


0.22 


0.25 


Z-score Ml duff/litter C 


-0.05 


2.7 


2.1 


2.2 


P-value Ml duff/litter C 


NA 


0.01 


0.05 


0.05 


Distribution duff/litter 


random 


clustered 


clustered 


clustered 


ANND tree 


NA 


0.61 


079 


0.59 


Z-score ANND tree 


NA 


3.7 


1.9 


2.9 


P-value ANND tree 


NA 


0.01 


0.05 


0.01 


Distribution tree 


NA 


clustered 


clustered 


clustered 



The spatial distribution of soil C pools is based Moran's / (Ml) for spatial 
autocorrelation. The spatial distribution is based on the Average Nearest 
Neighbor Distance (ANND) value. 



cover. These changes in vegetation cover lead to the loss 
of historical vegetation C; however, the loss of historical 
vegetation C is offset by large increases in carbon associ- 
ated with aboveground PJ growth and the accumulation 
of PJ surface-litter C. More notably, the changes in 
above ground carbon are highly variable across space in- 
dicating that the density of PJ woodlands in this region 
is the most important factor controlling the accumula- 
tion of C with woody encroachment. 



Woody encroachment and ecosystem C stocks 

Trees represent the major carbon stock in these SE Utah 
woodlands and the density and recent establishment of 
PJ woodland varied across our study area (Figures 2 
and 3). At the level of our experimental plots and in the 
larger analysis of aerial photographs, the spatial variation 
between relatively low and high density PJ woodland is 
responsible for the majority of the variation in C stock 
estimates within the Hop Creek study area. The net 
change in C between the historical condition and 
present states in the low-density and high-density plots 
respectively differs nearly threefold. This result high- 
lights the issues of field-based methods in sampling car- 
bon changes across this type of landscape. Had we 
simply sampled the LD and HD sites without attention 
to variation in density (e.g. in a long point-transect 
study), the mean value for these two sites would have 
been approximately 46% higher than our estimate for C 
based on analysis of orthophotos, which more accurately 
capture variation in tree density across this site. 
According to a one-correlation i-test power analysis, 258 
300 m plots delineated within Figure 3 would need to 
be sampled in order to sample 90% of the variation 
within the 25 ha area. 

The second major change in carbon stocks with PJ ex- 
pansion is in the litter/duff layer. The spatial distribution 
of litter/duff C follows the distribution of trees and is 
proportional to tree density (Tables 3 and 5). The accu- 
mulation of relatively recently fixed and non-mineral or 
aggregate stabilized C is similar to findings from a site in 
Southern Utah where the bulk of soil C accumulation 
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Table 6 Object-based image segmentation and ground 
based measurement comparisons for the LD and HD 
300 m 2 plots 



Measurement 


LD 


HD 


mage segmentation canopies 


32 


42 


Ground measured canopies 


J I 


A~) 


mage segmentation canopy 
area plot" 1 (m 2 ) 


137.00 


235.00 


Ground measured canopy area 
plot" 1 (m 2 ) 


1 30.05 


213.57 


mage segmentation mean 
canopy area tree" 1 (m 2 ) 


4.15 (0.55) 


5.46 (0.59) 


Ground measured mean 
canopy area (m 2 ) 


4.19 (0.69) 


5.08 (0.83) 


mage segmentation kg C 
plot" 1 


455.12 


846.1 1 


Ground based measured kg C 
plot" 1 using Huang et al. [49] 
calculation 


457.98 


833.39 


Ground based measured kg C 
plot" 1 using dbh 


327.26 


937.49 



Image segmentation tree C estimates per plot use the equation developed by 
Huang et al. [49]. Ground based tree C estimates are shown twice. Once using 
tree canopy area and the Huang et al. [49] equation and once using dbh and 
the appropriate tree species equation. 



following PJ establishment occurred in what is the 
equivalent to the duff /litter layer studied here [21]. Fur- 
ther, the estimated turnover time of this surface C pool 
in southern Utah was on the order of 10-20 years [21], 
a relatively rapid turnover time. The potential for rapid 
turnover of duff/litter combined with the vulnerability of 
these forests to drought-related mortality (as shown in 
the recent death of PJ stands in this study area) high- 
lights that the carbon gains associated with the expan- 
sion of these forests depend on the continued health and 
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Figure 4 Woody C increase area and magnitude from 1937 to 
2006 in the Hop Creek region. 



growth of these forests. Management activities that re- 
duce forest biomass or disturbances such as drought in- 
duced tree mortality or fire can lead to carbon losses 
that quickly reverse carbon accumulation in these 
forests. 

The carbon stock with the greatest inherent potential 
for long-term carbon stabilization is the mineral soil 
layer, where carbon could, in principle, be protected for 
decades to centuries. However, in this study, there were 
no appreciable differences in soil carbon stocks in areas 
with or without PJ forest. In these areas of SE Utah, the 
amount of carbon in the mineral soil layer is small. Add- 
itionally mineral layer soil C in soils similar to that of 
Hop Creek area that are mineral-stabilized only accounts 
for a small fraction of the soil C [21]. However, changes 
in soil carbon stocks with woody encroachment in other 
arid and semi-arid ecosystems differ from what was 
found in the Hop Creek area and suggests that soil C 
stock changes seem to be controlled by a number of 
confounding factors such as land-use history, past dis- 
turbance, and soil physical characteristics such as bulk 
density, texture, and mineralogy [17,19,34,35]. 

The temporal dynamics of carbon accumulation 

In the Hop Creek study area, the woodlands are still 
relatively young and appear to still be in a period of net 
C accumulation. Tree C accumulation rates for the plots 
examined here are rapid and increasing. For example, in 
the HD and LD woodland plots examined here, tree C 
accumulation averaged 0.37 kg tree" 1 year" 1 . The analysis 
based on image segmentation predicted an average tree 
C accumulation rate of 0.22 kg tree" 1 year" 1 ; however 
this estimate is based on the 69-year time-span derived 
from the 1937 and 2006 aerial photos. Using the average 
tree age from the plot-based data for the estimation 
period the rate of C accumulation for new trees using 
image segmentation is 0.30 kg yr" ; relatively similar to 
the plot-based estimate. At this rate of increase, these 
forests are accumulating C at an annual rate of 
~0.22 Mg ha" . It is notable that this study indicates that 
carbon accumulation in these forests is continuing with 
little sign of decline over time. This result differs from 
other studies that document a slowing of the rate of C 
accumulation as trees age and forests thicken [31,50,51]. 

The continuing accumulation of carbon in trees is 
striking when compared against evidence for widespread 
PJ mortality across the southwestern U.S. and in these 
experimental sites [32,33]. During the multi-year 
drought of 2000-2005, large areas of PJ woodlands died 
[32] including the site identified here as PJD. Just as 
woodland encroachment shifts the distribution of C 
within ecosystem compartments, woodland dieback also 
results in a redistribution of material. In this study, this 
shift is seen in the comparison of the LD and PJD plots, 
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which have similar duff/litter area cover but the PJD plot 
duff/litter layer is significantly deeper and richer in or- 
ganic C (Table 4). The timing of dieback (-2005) and 
the apparent difference in duff/litter C stocks of 0.11 kg 
C m~ 2 suggests that approximately 10% of tree C was 
transferred from live biomass to the duff/litter soil com- 
partment within a year of tree mortality. This observa- 
tion, combined with evidence for rapid (decadal) 
turnover rates of surface, un-stabilized C stocks in 
nearby PJ ecosystems [21] suggests that a sizeable frac- 
tion of the 20th century sequestration of C due to PJ en- 
croachment may be returned to the atmosphere due to 
recent tree-mortality in the next few decades. The com- 
bination of continuing carbon accumulation in growing 
PJ forests alongside patchy drought-related mortality fur- 
ther highlights the need to consider large-scale forest 
population dynamics in regional carbon stock estimates. 

Conclusions 

The importance of tree density and heterogeneous PJ 
drought mortality on C stocks further highlights the 
need for high-resolution, large spatial scale analyses of C 
mapping. However, the limitation of these techniques is 
that they cannot be used to assess understory plant 
cover and soil C stocks. Changes in surface duff carbon 
and soil carbon in arid land ecosystems tend to be highly 
variable and most easily quantified by plot scale mea- 
surements. Additionally, the loss of pre-encroachment 
vegetation cover (and carbon stocks) is difficult to estab- 
lish from aerial photos or other remote sensing analyses 
alone. In this site, if we relied solely on image segmenta- 
tion analysis, we would estimate that PJ encroachment 
since 1937 resulted in a woody C increase of 
15.36 Mg ha" 1 . However, our plot-based studies indicate 
that shrub C loss during woodland expansion accounts 
for approximately 2 Mg ha" of carbon loss thus lower- 
ing the carbon sequestration estimate by approximately 
13%. None-the-less, this study illustrates the importance 
of cross-scale high-resolution analyses of woodland bio- 
mass distribution and the sensitivity of ecosystem C 
stock estimates to variation in tree cover density. The 
distribution and density of woodlands is a critical deter- 
minant of landscape-scale C stocks and apparently more 
important, than variation due to different tree allometry 
equations. The counterpoint to continued C sequestra- 
tion in PJ ecosystems is the potential rapid and large loss 
of C due to large-scale mortality events. Our calculations 
suggest that as much as 25% of the C stabilized in the 
past century could be returned to the atmosphere in the 
next two decades in this regions woodland stands heavily 
impacted by drought-related mortality. Those losses 
could be larger if additional disturbances (e.g. wildfire) 
occur. 



Abbreviations 

ANN: Average nearest neighbor; BCA: Basal canopy area; C: Carbon; 
COFECHA: Computer-assisted tree-ring chronology composition system; 
ha: Hectare; MI: Moran's index; PJ: Pihon juniper; SOC: Soil organic carbon; 
RCD: Root collar diameter; N: Sample number. 

Competing interests 

The authors declare that they have no competing interest. 
Authors' contributions 

All authors were involved in developing the field work protocols and 
conducted field sampling. All authors read and approved the final 
manuscript. 

Acknowledgements 

This research was supported by the United States of Agriculture award 201 1- 
67003-30755 to J.N. and NASA North American Carbon Program Grant 
NACP-Asner-01. 

Author details 

1 Department of Geological Sciences and Department of Environmental 
Studies, University of Colorado, 2200 Colorado Avenue, Boulder, CO 80309, 
USA. department of Geography, National Taiwan University, Taipei 10617, 
Taiwan. 3 Department of Global Ecology, Carnegie Institution for Science, 260 
Panama Street, Stanford, CA 94305, USA. department of Ecology and 
Evolutionary Biology, University of Colorado, Boulder, CO 80309, USA. 

Received: 19 April 2013 Accepted: 22 July 2013 
Published: 5 September 2013 

References 

1. Scholes PJ, Archer SR: Tree-grass interactions in savannas. Annu Rev Ecol 
Syst 1997, 28:517-544. 

2. Archer S, Boutton TW, Hibbard KA: Trees in grasslands: biogeochemical 
consequences of woody plant expansion. In Global biogeochemical cycles 
in the climate system. Edited by Schulze ED, Heimann M, Harrison S, Holland 
E, Lloyd J, Prentice I, Schimel D. San Diego: Academic; 2001:115-138. 

3. Wessman CA, Archer S, Johnson LC, Asner GP: Woodland expansion in the 
US grasslands: assessing land-cover change and biogeochemical 
impacts. Land change science. In Observing, monitoring and understanding 
trajectories of change on the Earth's surface. Edited by Gutman G, Janetos AC, 
Justice CO. London: Kluwer Academic Publishers; 2004. 

4. Eldridge DJ, Bowker MA, Maestre FT, Roger E, Reynolds JF, Whitford WG: 
Impacts of shrub encroachment on ecosystem structure and 
functioning: towards a global sysnthesis. Ecol Lett 201 1, 14:709-722. 

5. Ratajczak Z, Nippert JB, Coliins SL: Woody encroachment decreases 
diversity across North American grasslands and savannas. Ecology 2012, 
93:697-703. 

6. Schlesinger WH, Reynolds GL, Cunningham LF, Huenneke WM, Jarrel RA, 
Whitford WG: Biological feedbacks in global desertification. Science 1990, 
247:1043-1048. 

7. Huenneke LF, Anderson JP, Remmenga M, Schlesinger WH: Deserification 
alters patterns of aboveground net primary production in Chihauhaun 
ecosystems. Glob Chang Biol 2002, 8:247-264. 

8. Hibbard KA, Schimel DS, Archer S, Ojima DS, Parton W: Grassland to 
woodland transitions: Integrating changes in landscape structure and 
biogeochemistry. Ecol Appl 2003, 13:91 1-926. 

9. Huxman TE, Wilcox BP, Breshears DD, Russel SL, Snyder KA, Small EE, Hultine 
K, Pockman WT, Jackson RB: Ecohydrological implications of woody plant 
encroachment. Ecology 2005, 86(2)308-31 9. 

10. Hughes FR, Archer SR, Asner AP, Wessman CA, McMurty C, Nelson J, Ansleyk 
JR: Changes in aboveground primary production and carbon and 
nitrogen pools accompanying woody plant encroachment in a 
temperate savanna. Glob Chang Biol 2006, 12:1733-1747. 

1 1 . Boutton TW, Liao JD, Filley TR, Archer SR: Beiowground carbon storage and 
dynamics accompanying woody plant encroachment in a subtropical 
savanna, Soil carbon and sequestration and the greenhouse effect. 2nd 
edition. Madison Wl, USA: Soil Science Society of America special 
publication; 2008. 57. 

12. Pacala SW, Birdsey R, Bridgham S, Caspersen J, Conant R, Davis K, Hales B, 
Houghton R, Jenkins J, Johnson M, Marland G, Paustian K, Socolow R, Tol R, 



Fernandez ef al. Carbon Balance and Management 2013, 8:8 
http://www.cbmjournal.eom/content/8/1/8 



Page 13 of 13 



Wofsy S: The North American carbon budget past and present, The first state 
of the carbon cycle report (SOCCR): the North American carbon budget and 
implications for the global carbon cycle. Departments of Energy and 
Commerce and NASA; 2007. 

1 3. Houghton RA, Hackler JL, Lawrence KT: The U.S. carbon budget: 
contributions from land-use change. Science 1999, 28(5427):574-578. 

14. Hunt GS, Moorcroft PR, Caspersen J, Shevliakova E, Houghton RA, Moore B 
II: Projecting the future of the U.S. carbon sink. Proc Nat Acad Sci United 
States 2002, 99(3):1 389-1 394. 

1 5. King AW, Hayes DJ, Huntzinger DN, West TO, Post WM: North American 
carbon dioxide sources and sinks: magnitude, attribution, and 
uncertainty. Frontiers Fco/2012, 10:512-519. 

16. Barger NN, Archer SR, Campbell J, Huang C, Morton JA, Knapp AK: Woody 
plant proliferation in North American srylands: a synthesis of impacts on 
ecosystem carbon balance. J Geophysical Res-Biogeosci 201 1, 1 16:G00K07. 

1 7. Jackson RB, Banner JL, Jobbagy EG, Pockman WT, Wall DH: Ecosystem 
carbon loss with plant invasion of grasslands. Nature 2002, 418:623-626. 

18. Knapp AK, Briggs JM, Collins SL, Archer SR, Bret-Harte MS, Ewers BE, Peters 
DP, Young DR, Shaver GR, Pendall E, Cleary MB: Shrub encroachment in 
North American grasslands: shifts in growth form dominance rapidly 
alters control of ecosystem carbon inputs. Glob Chang Biol 2008, 
14:615-623. 

1 9. Asner GP, Martin RE: Biogeochemistry of desertification and woody 
encroachment in grazing systems. In Ecosystems and land use change ft. 
Edited by Defries R, Asner GP, Houghton RA. Washington, DC: American 
Geophysical Union; 2004. 

20. Asner GP, Archer S: Livestock and the global carbon cycle. In Livestock in a 
changing landscape; drivers, consequences and responses. Edited by Steinfeld 
H, Mooney H, Schneider F, Neville LE. Washington, D.C: Island Press; 2010. 

21. Neff JC, Barger NN, Baisden WT, Fernandez DP, Asner GP: Soil carbon 
storage responses to expanding pihon-juniper populations in southern 
Utah. EcolAppI 2009, 19:1405-1416. 

22. Wheeler CW, Archer SR, Asner GP, McMurty CR: Climatic/edaphic controls 
on soil carbon/nitrogen response to shrub encroachment in desert 
grassland. Ecol Appl 2008, 17:1911-1928. 

23. Throop H, Archer R: Shrub (Prosopis velutina) encroachment in a semi- 
desert grassland: spatial-temporal changes in soil organic carbon and 
nitrogen pools. Glob Chang Biol 2008, 14:2420-2431. 

24. Liu F, Wu XB, Bai E, Boutton T, Archer SR: Quantifying soil organic carbon 
in complex landscapes. Glob Chang Biol 201 1, 17:1 1 19-1 129. 

25. Blackburn WH, Tueller PT: Pihon and juniper invasion in sagebrush 
communities in east-central Nevada. Ecology 1970, 51:841-848. 

26. Miller RF, Wigand PE: Holocene changes in semiarid pinon- 
juniperwoodlands. Bioscience 1994, 44:465-474. 

27. MILLER RF, ROSE JA: Fire history and western juniper encroachmentin 
sagebrush steppe. J Range Manage 1999, 52:550-559. 

28. Tausch RJ, Nowak RS: Fifty years of ecotone change between shrub and 
tree dominance in the Jack Springs Pihon research natural area. In 
Proceedings of the conference on Shrubland Ecotones. Edited by McArthur ED, 
Ostler WK, Wambolt CL. Ogden, UT, USA: US Department of Agriculture, 
Forest Service, Rocky Mountain Research Station; 1999. 

29. Wall TG, Miller RF, Svejcar TJ: Juniper encroachment into aspen in the 
Northwest Great Basin. J Range Manage 2001, 54:691-698. 

30. Weisberg PJ, Lingua E, Pillai RB: Spatial patterns of pihon-juniper 
woodland expansion in central Nevada. Rangeland Ecol Manag 2007, 
60:115-124. 

31. Shaw JD, Steed BE, DeBlander LT: Forest inventory and analysis (FIA) 
annual inventory answers the question: what is happening to Pihon- 
Juniper woodland? J For 2005, 6:280-2985. 

32. Breshears DD, Cobb NS, Rich PM, Price KP, Allen CD, Bailee RG, Romme WH, 
Kastens JH, Floyd ML, Belnap J, Anderson JJ, Myers OB, Meyer CW: Regional 
vegetation die-off in response to global-change type drought. 

Proc Natl Acad Sci USA 2005, 1 02: 1 5 1 44- 1 5 1 48. 

33. Shaw JD: Drought-related mortality in pihon-juniper woodlands: a test case for 
the FIA annual inventory system, Proceedings of the sixth annual forest 
inventory and analysis symposium; 2004 September 21-24; Denver, CO. 
Gen. Tech. Rep. WO-70. Washington, DC: U.S: Department of Agriculture 
Forest Service; 2005:65-71. 

34. Neff JC, Reynolds R, Belnap J, Lamothe P: Multi-decadal impacts of grazing 
on soil physical and biogeochemical properties in Southeast Utah. 
EcolAppI 2005, 15:87-95. 



35. Fernandez DP, Neff JC, Reynolds RL: Biogeochemical and ecological 
impacts of livestock grazing in semi-arid Southeastern Utah, USA, 
journal of Arid34 Environments. J Arid Environ 2007, 72:777-791. 

36. United States Department of Agriculture: Natural Resources Conservation 
Service. (ftp://ftp-fc.sc.egov.usda.gov/UT/Range/D35XY/035XY312UT.pdf). 

37. United States Department of Agriculture: Natural Resources Conservation 
Service, http://www.blm.g0v/pgdata/etc./medialib/blm/ut/richfield_f0/ 
planning/rmp/draft_rmp_eis/appendices.Par.29428.File.dat/Appendix%207. 
pdf. 

38. Huntoon PW, Billingsley GH, Breed WJ: Geologic map of Canyonlands 
National Park and vicinity. Utah: Canyonlands Natural History Assoc: 1982; 
scale 1:62500. 

39. U.S. Department of Agriculture Soil Conservation Service: Soil survey of 
Canyonlands area, Utah: parts of grand and San Juan Counties. Salt Lake City, 
Utah, USA: United States Department of Agriculture, Natural Resource 
Conservation Service; 1991. 

40. Darling ML: Structure and productivity of a pihon and juniper woodland in 
northern Arizona. Ph.D. disseration. Duke University, department of zoology; 
1967. 

41. Greir CC, Elliot KJ, McCullough DG: Biomass distribution and productivity 
of Pinus edulis-Juniperus monosperma woodlands of north-central 
Arizona. For Ecol Manag 1992, 50:331-350. 

42. Jenkins JC, Chojnacky DC, Heath LS, Birdsey RA: National-scale biomass 
estimators for United States tree species. For Sci 2003, 49:12-35. 

43. Sherrod LA, Dunn G, Peterson GA, Kolberg RL: Inorganic carbon analysis 
by modified pressure-calcimeter method. Soil Sci Soc Am J 2002, 
66:299-305. 

44. Stokes M, Smiley TL: Tree-ring dating. Chicago: The University of Chicago 
Press; 1968. 

45. Orvis KH, Grissino-Mayer HD: Standardizing the reporting of abrasive 
papers used to surface tree-ring samples. Tree-Ring Res 2002, 58:47-50. 

46. Holmes RL: Computer-assisted quality control in tree-ring dating and 
measurement. Tree-Ring Bull 1983, 43:69-78. 

47. Applequist MB: A simple pith locator for use with off-centre increment 
cores. J For 1 958, 56:141. 

48. Laliberte AS, Rango A, Havstad KM, Paris JF, Beck RF, McNelly R, Gonzalesv 
A: Object-oriented image analysis for mapping shrub encroachment 
from 1937 to 2003 in southern New Mexico. Remote Sens Environ 2004, 
93:198-210. 

49. Huang C, Asner GP, Martin RE, Barger NC, Neff JC: Multi-scale analysis of 
woody vegetation cover and aboveground biomass in pihon-juniper 
woodlands of the Colorado Plateau, USA. Ecol Appl 2009, 19:668-681. 

50. Lajtha K, Getz J: Photosynthesis and water-use efficiency in pinon-juniper 
communities along an elevation gradient in northern New Mexico. 
Oecologia 1993, 94:95-101. 

51. Ryan MG, Blinkley D, Fownes DH: Age related decline in forest 
productivity: pattern and process Adv. Ecol Res 1997, 27:213-262. 



doi:1 0.1 186/1 750-0680-8-8 

Cite this article as: Fernandez ef al:. Twentieth century carbon stock 
changes related to Pinon-Juniper expansion into a black sagebrush 
community. Carbon Balance and Management 201 3 8:8. 



Submit your next manuscript to BioMed Central 
and take full advantage of: 

• Convenient online submission 

• Thorough peer review 

• No space constraints or color figure charges 

• Immediate publication on acceptance 

• Inclusion in PubMed, CAS, Scopus and Google Scholar 

• Research which is freely available for redistribution 



Submit your manuscript at 
www.biomedcentral.com/submit 



o 



BioMed Central 



